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LANE FORMATION BY SIDE-STEPPING 


MARTIN BURGER*, SABINE HITTMEIRt, HELENE RANETBAUER*, 
AND MARIE-THERESE WOLFRAM* 


Abstract. In this paper we study a system of nonlinear partial differential 
equations, which describes the evolution of two pedestrian groups moving in 
opposite direction. The pedestrian dynamics are driven by aversion and cohe¬ 
sion, i.e. the tendency to follow individuals from the own group and step aside 
in the case of contraflow. We start with a 2D lattice based approach, in which 
the transition rates reflect the described dynamics, and derive the correspond¬ 
ing PDE system by formally passing to the limit in the spatial and temporal 
discretization. We discuss the existence of special stationary solutions, which 
correspond to the formation of directional lanes and prove existence of global 
in time bounded weak solutions. The proof is based on an approximation ar¬ 
gument and entropy inequalities. Furthermore we illustrate the behavior of 
the system with numerical simulations. 


Keywords: diffusion, size exclusion, cross-diffusion, global existence of solutions. 

1. Introduction 

In the last decades demographics, urbanization and changes in our society re¬ 
sulted in an increased emergence of large pedestrian crowds, for example the com¬ 
muter traffic in urban underground stations, political demonstrations or the evac¬ 
uation of large buildings. Understanding the dynamics of these crowds has become 
a fast growing and important field of research. The first research activities started 
in the field of transportation research, physics and social sciences, but the ongoing 
development of mathematical models initiated a lot of research also in the applied 
mathematics community. Nowadays mathematical tools to analyze and investigate 
the derived models provide useful new insights into the dynamics of pedestrian 
crowds. 

A variety of different mathematical models has been proposed in the past which 
can be generally classified into microscopic and macroscopic approaches. In the 
microscopic framework the dynamics of each individual is modeled taking into ac¬ 
count social interactions with all others as well as interactions with the physical 
surrounding. This approach results in high dimensional and very complex systems 
of equations. Examples include the social force model by Helbing (cf. [5], [T7], [TB]1. 
cellular automata (cf. [12], |3], [11], [T]) or stochastic optimal control approaches, 
cf. pH] . 

Macroscopic models, where the crowd is treated as a density, can be derived by 
coarse graining procedures from microscopic equations (see e.g. [6]), leading to 
nonlinear conservation laws or coupled systems of such (see e.g. [H], [TU], [5]). 
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Other approaches heuristically motivating macroscopic models are based upon op¬ 
timal transportation theory, cf. [5^, mean field games (cf. [23], [5]) or optimal 
control cf. [14]. Piccoli and co-workers (cf. [28] and [11]) proposed a measure based 
approach capable to describe pedestrian dynamics on both microscopic and macro¬ 
scopic scale - hence bridging the gap between the two description levels. Recently 
there has been an increasing interest in kinetic models and their respective hydro- 
dynamic limits in pedestrian dynamics, see for example m and m- 
For an extensive review on the mathematical literature concerning crowd dynamics 
and the closely related field of traffic dynamics we refer to [2] . 

In this paper we (formally) derive and rigorously analyze a PDF system describing 
the evolution of two pedestrian groups moving in opposite direction. The indi¬ 
vidual dynamics are driven by two forces, cohesion and aversion. We show that 
this minimal dynamics already result in complex macroscopic features, namely the 
formation of directional lanes. We start with a 2D lattice model, in which the tran¬ 
sition rates, i.e. the rate at which a particle jumps from one site to the next, express 
the tendency of individuals to stay within their own group (i.e. follow individuals 
moving in the same direction) while stepping aside when individuals from the other 
group approach. The corresponding mean-field PDF model can be derived by a 
Taylor expansion (up to second order) and is a nonlinear cross diffusion system 
with degenerate mobilities. 

Similar models have been proposed in the literature, for example in the context of 
ion transport, cf. [3] or population dynamics, cf. [HD]. The coherent difference of our 
model to these works are additional challenging features, namely a perturbed gra¬ 
dient flow structure as well as an anisotropic degenerate diffusion matrix. Although 
the system lacks the classical gradient flow structure, we can show that the entropy 
grows at most linearly in time. The corresponding entropy estimates are a crucial 
ingredient for deriving the global existence result for bounded weak solutions. The 
existence proof is based on an implicit time discretization and an 77^-regularization 
of the time-discrete problem. Note that we follow a different approach than Jiingel 
in [21], which has the advantage that the method is based on an 77^-regularization 
only and does not require the additional Bilaplace operator. We define a fixed point 
operator in L'^{Q) and use Schauder’s fixed point theorem to deduce the existence 
of a solution to the regularized problem. The derived entropy estimates as well as a 
generalized version of the Aubin-Lions lemma justify the limit in the regularization 
parameter. 

This paper is organized as follows: In Section [2] we present the 2D lattice based 
model and derive its (formal) mean-field limit via Taylor expansion up to second 
order. Furthermore we discuss the existence of special stationary solutions in Sec¬ 
tion [Ql Section [3] focuses on structural features of the resulting PDF system, 
such as the corresponding entropy functional and the related dissipation inequality. 
Furthermore we study the boundedness of the densities, which is an essential pre¬ 
requisite for the global existence proof outlined in Section |3J Finally we illustrate 
the behavior of the model with various numerical experiments, which reproduce 
well known phenomena such as lane formation in Section |SJ 

2. Mathematical modelling 

In this section we present the formal derivation of the proposed PDF model 
from a microscopic discrete lattice approach. We consider two groups of individuals 
moving in opposite direction, i.e. one group is moving to the right, the other to the 
left. The individual dynamics are driven by two basic objectives: first individuals 
try to stay within or close to their own group, i.e. pedestrians walking in the same 
direction. Moreover they step aside when being approached by an individual moving 
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in the opposite direction. Based on this minimal interaction rules we derive the 
corresponding PDE model by Taylor expansion up to second order in the following. 


2.1. The microscopic model. Throughout this paper we refer to the groups of 
individuals moving to the right and left as red and blue individuals respectively. 
Their dynamics are driven by the objectives described above and correspond to 
cohesion and aversion. Let us consider a domain C partitioned into an 
equidistant grid of mesh size h. Each grid point {xi,yj) = {ih,jh), i £ {0,.. .N} 
and j £ can be occupied by either a red or a blue individual. The 

probability to find a red individual at time t at location {xi,yj) is given by: 

rij{t) = P(red individual is at position (xi,yj) at time t), 

with an analogous definition for We set tk = kh tor k £ N and use the abbre¬ 

viation Xij = rij{tk) if the time step tk is obvious. The dynamics of the individuals 
are driven by the evolution of the probabilities r and b. These probabilities depend 
on the transition rates of individuals. Let ol denote the transition rate 

of an individual to move from the discrete point {xi,yj) to (xi+i,yj). We define 
the transition probabilities for the reds as: 

= (1 _ + ar,+2.,), 

^ ( 1 ) 


where p = r -b 6 , 0 < 70 , 71,72 < 1 and 0 < a < ^. The factor (1 — p) corresponds 
to size exclusion, i.e. an individual cannot jump into the neighboring cell if it is 
occupied. Note that we assume that individuals only anticipate the dynamics in 
their direction of movement, i.e. they do not look backwards, which is reasonable 
when modeling the movement of pedestrians. The second factor in the transition 
probabilities m corresponds to cohesion and aversion. If a > 0 the probability of 
moving in the walking direction is increased if the individual in front, i.e. at position 
(xi+ 2 ,?/j), is moving in the same direction (assuming that the cell {xi+i,yj) is not 
occupied). 

Aversion corresponds to sidestepping. If 71 > 72 > 0, an individual steps aside if 
another individual, in m a blue particle located at {xi+i,yj), is approaching. If 
7i > 72 , there is a preference to make a step to the right hand side with respect 
to their direction of movement, if 72 > 71 , to the left. From the perspective of an 
observer red individuals prefer to make a jump down if a blue individual is ahead 
of them in the case 71 > 72 . The parameter 70 > 0 includes diffusion in the y- 
direction. In the case of no diffusion, i.e. 70 = 0, individuals only step aside when 
being approached by an individual moving in opposite direction. 

The master equation for the red particles then reads as 


- + 7-fe}^bJ+i}^ n^Atu). 


The probability to find a red particle at location (xi,yj) in space corresponds to 
the probability that a particle located at (xi-i^yj) jumps forward (first term), 
particles located above or below, i.e. at {xi,yj±i) jump up or down (second line), 
minus the probability that a particle located at (xi, yj) moves forward or steps aside 
(third line). The corresponding transition rates for the blue particles are defined 
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accordingly to m by: 


— (1 Pi—“1“ ^ 1 


'h 

= (1 _ 


% 




= (1 - Pi.i-l)(70 + 72 n-i,j)- 


The master equation for the blue particles has the same structure as i.e.: 
bij{tk+i) = bij{tk) + 


+ '6 Oi^j^l[tk) 

( 


-it; 


b 

{tjI-s-It-Ij} 




% 




bij{tk)- 


( 3 ) 


(4) 


2.2. Derivation of the macroscopic model. In the following we shall consider 
the formal limit h — At = Ax = Ay —>■ 0 in equations ([U and (|3]) to derive the 
corresponding PDE system. After performing a Taylor expansion up to second 
order, we obtain 

dtr = -V ■ Jr, 
dtb = -V • Jb, 

where 

1^(1 - p)(l + ar)r + | [aa;(r(l - p)(l + ar)) - 2((1 - p)da;r)]\ 
-(7i-72)(l-p)6r-|[(7i+72)((l-p)5y(r6)+6rdyp) ’ 

V +270 ((l-p)dj/r + r5yp) +2(71-72)(l-p)r5,,&] J 

and 

''-{1 - p)(l + ab)b + I [d,r{b{l - p){l + ab)) - 2((1 - p)9a;6)]\ 
(7i-72)(l-p)6r-|[(7i+72)((l-p)ay(r6)+6raj^p) 

\ +270 ((1 - p)dyb + bdyp) + 2(71 - 72)(1 - p)bda,r] ) 

denote the fluxes for r and b respectively. The first order terms correspond to the 
movement of the reds and blues to the right and left in x-direction respectively as 
well as to the preference of either stepping to the right or left in y-direction (depend¬ 
ing on the difference 71 — 72 )- The second order terms correspond to the cross diffu¬ 
sion terms where the prefactor h denotes the lattice size. We consider system on 
n X (0, T), where C is a bounded domain. In our computational examples, see 
SectionjHl the domain 17 corresponds to a corridor, i.e. 17 = [—La,, La,] x [— 
with Ly Lx- As individuals cannot penetrate the walls, we set no flux boundary 
conditions on the top and bottom, i.e. 





= 0 


at ?/ = ALy . 


At the entrance and exit of the corridor, i.e. at x = ALx, we assume periodic 
boundary conditions. Note that Robin type boundary conditions, where the in- and 
outfluxes at the entrance and exits are directly proportional to the local density, 
would be more realistic. The boundary conditions set above correspond to the 
simplest choice and shall serve as a starting point for the investigation of more 
realistic and complex models in the near future, cf. [7]. 

We would like to remark that the lengthy Taylor expansion and formal limiting 
procedure can be accomplished automatically using computer algebra techniques, 
even for more general classes of models, see P5] . 
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2.3. Stationary Solutions. In this last part of the modelling section we study the 
existence of specific stationary solutions, which correspond to the formation of lanes. 
These segregation phenomena can be observed in crowded streets with pedestrians 
as well as in experiments. Lane formation is a rather intuitive phenomenon, but a 
strict mathematical definition is less obvious. In the following we shall distinguish 
between strict segregation and the case when still some pedestrians might get into 
the counterflow, leading to the definition: 


Definition 1. Let {r,b) denote a stationary solution to system @ for 71 > 72 , 
which is x-independent, i.e. for all x,xq € [—Lx^Lx] and any y € [—Ly,Ly] we 
have (r,b)(x,y) = {r,b){xQ,y). Considering therefore {r,b) as a function of y only, 
we call {r,b) S L°°[—Ly,Ly] x L°°[—Ly,Ly] 


• a solution with strong lane formation, if the functions r and b have a com¬ 
pact support in y-direction with 

supp(r) n supp( 6 ) = 0 and sup {supp(r)} < inf {supp( 6 )}. 

y&[-Ly,Ly\ y(^[-Ly,Ly] 

• a solution with weak lane formation, if the sufficiently smooth solution (r, b) 
satisfies 

dyr < 0 and dyb > 0 , 
and there exists a point y G {—Ly, Ly), such that 

r{y) = b{y). 


Note that the definition of weak lane formation has to be changed accordingly if 
individuals have the preference to step to the left instead of right, i.e. 72 > 7i- 
We expect that the side-stepping initiates the formation of directional lanes, an 
assumption that has also been confirmed by the numerical experiments in Section 
[ 5 ] for specific ranges of parameters. In particular we consider system ([Sl) in the case 
0 < a < i and 71 — 72 = 0{h), i.e. 

dtr = -dx ((1 - p)(l -b ar)r) +hdy ((1 - p)br) - ^dl{r{l - p)(l -b ar)) 

-b hdx{{l - p)dxr) -b ^(71 + l 2 )dy ((1 - p)dy{rb) + brdyp) 

+ h')Qdy ((1 - p)dyr + rdyp) 

h 

dtb = dx {{1 - p){l + ab)b) - hdy ((1 - p)br) - -dl{h{l - p){l -b ab)) 

-b hdx{{l - p)dxb) + ^(71 -b "i2)dy ((1 - p)dy{rb) + hrdyp) 

-b h'jody ((1 - p)dyb -b bdyp ), 


where we set without loss of generality 71 —72 = h. Note that the second order terms 
in h are dropped out, i.e. the terms including x- and y-derivatives are neglected. In 
this case we can prove weak lane formation for 70 > 0 and postulate the formation 
of strong lanes as 70 —>■ 0 . 

Therefore, we consider system ([5]) and analyze its equilibrium solutions which are 
constant in x-direction. In this case system (jS]) reduces to 


0 = (1 - p)rb + ((1 _ p)dy{rb) + brdyp) + 70 ((I - p)dyr + rdyp ), 

0 = -(1 - p)rb + ((1 _ p)dy{rb) -b brdyp) -b 70 ((1 - p)dyb + bdyp). 


(7a) 

(7b) 


2 
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Note that we have assumed a preference for stepping to the right in ([5]), which 
corresponds to the different sign in the first terms of ©• If p < 1, we can rewrite 
(HI) as 


0 = 


0 = -: 


rb 

1- p 
rb 


7i +72 


dy 


7i + 72 


1-p 


dy 


rb 

1-p 

rb 


■7o9, 


r 


1-p 


lod. 


" ' l-p 
b 


l-p 


Summation of (I5al) and (IFbl) and subsequent integration gives 

/ \ 'I'b P ^ 

7i +72 -j-h7o:j-= C, 

l-p l-p 


(8a) 

(8b) 


(9) 


for some constant C. Equation (|S]) allows us to study the behavior of stationary 
solution curves with respect to the densities r and b. Figure [T] illustrates these 
stationary solutions in the case 71 — 72 = I* for different values of C. If r = 0 



Figure 1. Stationary solution curves for 71 = 0.5, 72 = 0.4 and 
70 = 0.001 

or 5 = 0 then dyr = 0 01 dyb = 0 respectively. Hence solution curves can get 
arbitrarily close to the r- and 6 -axes, but they can only reach them in the case of 
a trivial solution curve, i.e. consisting only of one stationary point lying on one of 
the axes. The actual starting and end points of the solution curves as well as the 
corresponding constants C depend the chosen parameters and on the initial masses 
of the system, i.e. on 

My := / rdxdy and Mf, := / bdxdy. 

Jq Jn 

In the case of small values of C we observe a quick change of the densities r and 6 
from high to low values and the other way around. For larger values the densities 
increase or respectively decrease slower along the solution curves. 

The following additional solution properties can be deduced from equations ([5]) and 

(Ell- 

Lemma 1. Let {r,b) denote solutions to system ([8]) and let C S R'*" be a constant 
with 0 < (7 < 1. 

(i) There exists no solution {r,b) with p = C and r,b > 0. 

(ii) There exists no solution {r,b) with r = 6 > 0. 
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(in) Any solution (r, b) is monotone with dyU < 0 and dyb > 0. 


Proof. To show Q we assume to the contrary that there exists a solution with 
p = C and r,b > 0. Then equation m implies that rb is a positive constant and 
therefore the same holds true for r and b individually. This is a contradiction to 
(0) as is only true if rb = 0 . 

To see ([u]) we again argue by contradiction. If r = 6 we immediately deduce 
from system (|5]) that r = b = 0. 

To prove the monotinicity properties in (Hull we first observe that (I 8 all and (I 8 bll 
imply 


+ Cir'^ 


< 0 


and 


(^) 


> 0 , 


( 10 ) 


for some constant Ci > 0. This allows to exclude the existence of a y G [~Ly, Ly] 
with dyr{y) = 0 , since in this case equations (ITUl) would yield dyb{y) < 0 as well 
as dyb{y) > 0. Therefore r has to be monotone and due to symmetry b is also 
monotone with the opposite sign. 

To show the stated signs of the derivatives we assume dyU > 0 and dyb < 0. 
Subtracting the equations in (nni then leads to {r — b)dyp < 0 and thus to a 
contradiction, since for dyp > 0 equation ( 1101 ) as well as the assumptions imply 
r — b > 0, whereas for dyp < 0 the second equation in (ITUl) gives r — b < 0. We 
therefore obtain the desired monotonicity properties dyr < 0 and dyb > 0 . □ 


Lemma [T] indicates the existence of weak lane formation. From (Iml) we know 
that r and b are monotone functions which are strictly positive. Hence there exists 
a single point y € {—Ly,Ly) where r{y) = b{y). Due to the side-stepping tendency 
the reds will move to the bottom, while the blues move up. In the case of equal 
masses it is impossible that one density is larger than the other on the whole domain, 
which implies the formation of weak lanes in the sense of Definition [1] in this case: 

Theorem 1. Let 70 > 0, 71 > 72 and Mr = M}, = M. Then system ([ 8 ]) has non¬ 
trivial stationary states, and any stationary solution constant in the x-direction 
exhibits weak lane formation. 

Further properties of solutions to © and m can be observed for different as¬ 
ymptotic parameter regimes: 

• 7 o —>■ 0: We deduce from m and m that rb —> 0 , i.e the smaller 70 , the 
sharper the separation of r and b. 

• 7 i + 72 —>■ 0: In this case p = (7 for some constant 0 < C < 1 and 
dyr = Cir{C — r) for some constant Ci > 0 which corresponds to lane 
formation. 


Remark 2. If pedestrians have the preference to step to the left instead of to the 
right the monotonicity behavior of r and b is reversed. 

3. Basic properties 

In this section we discuss basic properties of system (0). In the following we set 
7o > 0, 7 := 7i = 72 and a = 0. Then system © reads as 

dtr = -d,r ((1 - p)r) + ^d^ ((1 - p)da;r rd^^p) 

+ h [lody ((1 - p)dyr + rdyp) P jdy ((1 - p)dy{rb) + brdyp)] 

dtb = dx ((1 - p)b) -k ((1 - p)d,rb -k bd^p) 

+ h [lody ((1 - p)dyh + bdyp) -k 79 j, ((1 - p)dy{rb) + brdyp)]. 


( 11 ) 
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We shall prove global existence of weak solutions of system (HH) in Section H a 
result which can be extended to the case a > 0 and 71—72 = 0{h). The proof 
uses several structural features of system m, such as the corresponding entropy 
functional and the boundedness of solutions, which we discuss in this section. 


3.1. Entropy functional. A key point in the existence analysis are estimates 
based on the corresponding entropy functional 


£ := 


/ r(\ogr — 1) + b(logb — 1) dx dy 

Jn 

+ ^ ^(1 - /3)(log(l - p) - 1) + ^rVr + ^bVbdxdy, 


( 12 ) 


where the potentials Vr{x,y) = —x and Vb(x,y) = x correspond to the motion of 
the red and blue individuals to the right and left respectively. Note the difference 
in the prefactor ^ of the entropy term (1 — p)(log(l — p) — 1) compared to other 
entropies used in the literature for similar PDE models, cf. a, m- This prefactor 
results from the anisotropic diffusion as we shall explain in the following. 
Introducing the entropy variables u and v 

12 12 

U := drS = logr- - log(l - p) + -K and v := dhS = \ogh-- log(l - p) + -14, 


h 


we can rewrite dm as follows 



/ 



/ ia.p\\ 

{dtr\ h (V 0\ 
\dtb) 2 Vo V' 

M 

dyU 

dxV 

+ 

lordyp 


\ 

\dyv) 


\lobdypJ ) 


where 


M = 


/(l-p)r 
0 
0 
0 


0 


0 


(13) 


\ 


270(1 — p)r + 27(1 — p)r& 0 27(1 — p)r& 

0 (l-p)& 0 

27(1 — p)r6 0 270(1 — p)5 + 27(1 — p)r6y 

We observe from equation m that we do not have a gradient flow structure. 
The additional terms result from the different structure of the second order terms 
in (HH). They are either of the form r(l — p), 6(1 — p) or r6(l — p), which correspond 
to different entropies. This lack of structure results in the different prefactor in 

(na. 

Note that the entropy functional (na is also not an entropy in the classical sense as 
we cannot ensure that it is non-increasing. Nevertheless the entropy grows at most 
linearly in time, which is sufficient for proving existence of global weak solutions. 


Lemma 2. Let r, 6 : 17 — ^ 6 e o sufficiently smooth solution to system ca 

satisfying 

0 < r,b and p < 1. 

Then there exists a constant C > 0 such that 

^+Vo<C, (14) 

where 

■Do = Cq [ (l-p)|Vv^p + (l-p)|vy6p + |Vv^Wp + |Vppdxdp, 

Jn 

for some constant Co > 0. 
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Proof. System (fT^ enables us to deduce the entropy dissipation relation: 

df 


dt 


= / {u dtr + V dtb)dx dy 

Jo. 




Vu\ fVu 
Vv) Ivu 


( \ 

lordyp 

^doop 

V lodyp J 


Vu 

Vu 


dx dy 


h 


(15) 


(1 - p){r{dxuf + b{dxvf) + -d^p{rdxU + bd^v) 

+ 270 - P){'fi.dyuf + b{dyvY) + ^dyp{rdyU + bdyv) 


—p)rb{{dyU +dyvY dxdy, 

where we have used integration by parts. For the non-quadratic term in x-direction, 
we use the fact that 


d^r dxP 

OxU = -h 

r 


_2 _ dopb ^ dxP ,2 

2(1 — p) h' b 2(1 — p) h' 


and deduce that 
1 

~2 


J dxpirdxU + bd^v) dxdy = -^ J ^1 -b 2(^1^- p) ) {9xp)'^dxdy 

+ T I 9xp{r - b) dx dy. 

“ Ja 


(16) 


The first term on the right hand side is negative, for the second we derive that 
1 


■ f d^p{r -b)dxdy=\ / (1 - p){dxr - d^b) 
Jn P JQ 


rdxU - bdxV - ^ 


= y (1 - p){rdxU - bdxv) ^ J 9 xp{r - b) dx dy 
J p{l- p) dxdy. 


2 

7)2 


Therefore we obtain 


^ J dxp{r-b)dxdy = ^ J {1 - p){rdxU - bd^v) dx dy +-^ J p(l - p) dxdy. 

As 0 < p < 1 and the integration is over a bounded domain, there is a positive 
constant C such that 

P(1 - P) dxdy < C. 

Applying Young’s inequality, we get 

f {1 — p)rdxU dx dy < — f {l — p)rdxdy + — [ {1 — p)r{dxuf‘dx dy 
Jn b Jq 4 Jq 

and therefore 

I {I - p){rdxU - bdxv) dx dy < ^ [ p{l-p)dxdy 

Jn b Jq 

+ 7 f {I- p){r{dxuf‘ + b{dxvf) dxdy. 

4 Jn 
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Altogether we deduce the following estimate from (ITBl) : 

I da:p{rd^u + bdxv)dxdy <-]- / [l+ ^ 


{d^pf dx dy 


2(1-P), 

+ ^ [ {1-p){r{dxu)'^+ b(dxvf)dxdy 
6 Jn 

+ ^ P)dxdy. 


We use the same arguments for the term — 5 /q dyp{rdyU + bdyv) dx dy and obtain 
the following entropy dissipation from (1151) : 


f 

— = / {udtT + V dtb)dx dy 

dt Jn 


< — 


^(1 - p){r{d^uf + b{d^vf) + ^ 


2(1-P) 


{{d^pY + 2-fo{dypf) 


- (1 + 27 o)^p (1 - p) + ^7o(l - p){r{dyuY + b{dyvf ) 


+ 27(1 - p)rb{{dyU + dyvY 


dx dy < C. 


For the analysis it will be sufficient to use a reduced version of the entropy inequality, 
given by 

^ < -Co / (1 - p)(r|Vup + 5|VuH + IVpp dxdy + C =: Pq, 

Jn 

where Co := ^min(l, 7 o). Using the definitions of u and v, applying Young’s 
inequality to estimate the mixed terms as well as the fact that 

2 


r(l - p) 


V log 


2(1-P) 


+ 5(l-p) 


V log 


= 4(l-p)|VvC| +4(l-p) Vv^ +4p Va/1-p +2|Vp| 


2(1-P) 
2 


we obtain 


'Do<-Co [ 2(1 - p)|VvCp + 2(1 - p)|vV5|2 + 2p|vyi - pp + |Vpp dxdy 

Jn 

+ 4^0 f Y-p){r\^Vr\‘^+b\VVb\‘^) dxdy-Co f \V pl"^ dx dy + C. 

Since |VUrp = |VVhp = 1 and 2p\V^Jl — p\^ + |Vpp > 2|Vv'l — pP, we get the 
estimate 

< -Co / (1 - p)(|VvCp + IVv^n + ivywp + |Vpp dxdy + C (17) 
Jn 

for some constant C > 0, which concludes the proof. □ 


dt 


3.2. Positivity. We want the global weak solution of system m to satisfy 0 < 
r{t), b{t), p{t) < 1 for all t > 0 , if the latter condition is prescribed for the initial 
data. System (ED can be written in the form 






/(i-p)rvv;\\ 

\{l-p)b^Vb)) 
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where A = A{r, b) is the diffusion matrix given by 

— b) 0 r ^ \ 

0 27(1 — 6)6 + 270(1 — 6 ) 0 27(1 — r)r + 27 or 

6 0 (1 -r) 0 

0 27(1 — 6)6 + 2706 0 27(1 — r)r + 270(1 — r)y 

Note that the diffusion matrix is neither symmetric nor positive definite in gen¬ 
eral. Hence, we cannot use the maximum principle to prove nonnegativity and 
boundedness of r, 6 and p. However, the system allows to use a more direct ap¬ 
proach to deduce upper and lower bounds for the variables r, 6 and p. We therefore 
consider the entropy density 

^ R, r(logr-l)+6(log6-l) + i(l-p)(log(l-p)-l) + ^rV; + ^&f4, 

where 

M = :r>0,b>0,r + b<lY (18) 


/(I 


A{r, 6) = 


Lemma 3. The function E : A4 is strictly convex and belongs to C^{A4). Its 
gradient DE : At —>■ zs invertible and the inverse of the Hessian D^E : A4 —> 

is uniformly bounded. 

Proof. The invertibility of DE can be shown directly. Using the definitions of the 
entropy variables u and v, we get 

T 4x vh 

u — V = log --— and u + v = log-. 

b h r — p 


Solving these relations for r gives 


r = 6e '> 


and r = 


(1 - 6)e“+“ 


6 + e“+" ’ 

which leads to a quadratic equation in 6 with exactly one positive solution 

, 2 


6 = 6(u,u) =-i (e“+’' + e-Te2’') 


V 






and therefore 


(e“+" + eTe2“) + 

/^fe“+’' + eTe2“V \ 

--- — + €"^6^“ 


1 ) 


Simple calculations ensure that G M. 

To show the uniform boundedness of the inverse of D^E : Ai 
that 


we observe 


DE 


and D^E ( ^ 
6 


drU dbU 
drV dbV 


^ r ^ 2(l-p) 

\ 2(1^ 


1 p-"! 

U “T 


b ^ 2(l-p)^ 


Since 0 < r, 6, p < 1, we can deduce that the inverse of D‘^E exists and is bounded 
in At. □ 
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Hence, LemmalHlensures that if there exists a weak solution (w, v) G T; 


to (EH, the original variables ( , ) = (DE) 




satisfy 




G Af for t > 0 


almost everywhere. This gives us L°°-bounds, necessary for the global in time ex¬ 
istence proof in the following section. 


4. Main result 

We start this section by stating the notion of weak solutions to system EH). 


Definition 2. A function (r^b) : x (0,T) ^ M is called a weak solution to 

system EH if it satisfies the formulation 





70 


da;r{l - p) + rda;p\ f 
doob{l - p) + bdxp) 

dyr{l - p)+ rdyp\ _ 

dyb{l - p) + bdypJ 


dy<^>l 


dy{rb){l-p) 
dy{rb){l- p) 

(l-p)rVW\ 
(i-p)6vyJ ■ 


. rbdyp'\ 

.rbdypj 


dy^i 

dy^2 


dx dy dt 
dxdydt (19) 
dx dy dt 


dx dy dt = 0, 


for all ^ 1,^2 G 


Theorem 3. (Global existence) Let T > 0, let {ro,bo) : LI —^ A4, where A4 is 
defined by EH, be a measurable function such that E{ro,bo) G L^{Ll). Then there 
exists a weak solution {r,b) : fl x (0,T) —>■ At in the sense of (1191) with periodic 
boundary conditions in x-direction and no-flux boundary conditions in y-direction 
satisfying 

dtr, dtbG L'^{0,T;H\ny), 
p, €L^{0,T-,H\n)), 

{l-p)W^/f, {l-p)wVb G L^{ 0 ,T;L'^{n)). 

Moreover, the weak solution satisfies the following entropy dissipation inequality: 

^+Pi<C, (20) 

dt 

where 

Vi = Co [ (l-p)2|Vv^p + (l-p)2|vC6|2-|-|Vv^Wp + |Vppdxdy 

JQ 

and Co and C are the constants from 

We would like to mention the different dissipation term in (l20ll . In particular, 
since the convergence properties are not strong enough to pass to the limit in the 
entropy dissipation El, we obtain a modified entropy inequality (j20p . 

A major difference in the analysis of the system compared to related ones in the 
literature (cf. [30], [H, EH, EH) is the fact that we have an anisotropic diffusion and 
no gradient flow structure, which requires a different entropy and a priori estimates. 

The following existence proof is based on an approximation of EH- The basis 
of the approximation argument is the following formulation of system EH: 
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( 21 ) 


where 



/(I - p)r(l + ^r) 

0 

(l — p)rb 

2-p 

0 \ 

h 

0 

2(l-p)r(jo(§E^)+7b) 

0 

2(1 - p)r 6 (^ 4 - 7 ) 

2 

(1 - p)b(l + 

0 

(l — p)rb 

2-p 

0 


\ 0 

2(1 - p)b(7o(§E^) + jr) 

0 

2(l-p)&(^+7r)/ 


and 


H = 


0 

V 0 


The positive semi-definiteness of the matrix G(r, b) can be proven using a similar 
approach as we have seen in Subsection 13.11 

We discretize system m in time using the implicit Euler scheme with time step 
r > 0 which results in a recursive sequence of elliptic problems. These are modified 
by adding higher order regularization terms. The corresponding weak formulation 
is given by 


f^k- rk-i\ 
T \ bk - bk-i j 


Auk + 

\Avk +Vk J ' 


( 22 ) 


The regularization guarantees coercivity of the elliptic system in This is 

in contrast to |30j who used a stronger regularization by introducing a Bilaplacian. 
The existence proof is divided into several steps. First we show existence of weak 
solutions to the regularized, discrete in time problem by applying Lax-Milgram to 
a linearized version of the problem (1221) and using the Schauder fixed point theorem 
to conclude the existence result for the corresponding nonlinear problem. 

Finally uniform a priori estimates in r and the use of a generalized Aubin-Lions 
lemma (cf. [3D]) allow to pass to the limit r —^ 0. Note that one can also use the 
Kolmogorov-Riesz theorem in a similar fashion to [Ij. 


4.1. Time discretization and regularization of system dTT|) . We start by 
studying the regularized time discrete system. Recall that the entropy variables 
are defined as {u,v) = DE{r,b) for (r, 6) G A4. Lemma |3| ensures that DE is 
invertible, hence we set (r, 6) = {DE)~^{u,v) for {u,v) G 

Let T > 0, N G N and let r = T /N be the time step size. We split the time interval 
into the subintervals 


rp 

Then for given functions (r^-i, bk-i) G M, which approximate (r, b) at time T(fc — 
1), we want to find {rk,bk) G M solving the regularized time discrete problem (1331) 
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in the weak formulation: 


T Jn - bk-i) 



dx dy + 


+ / H{rk,bk) 

JQ 


/V$i 

\y<i>2 



T 

G{rk,bk) 


/'Vuk\ 

\yvk j 


dx dy 


^ dxdy + tR 




for ($i,$2) G H^{Vt) X where {rk,bk) = DE ^{uk,Vk) and 


( 23 ) 


(($0 ’ (^fc)) ^ J ^1^^= + + V$i • Vufc + V$2 • Vufc dxdy. 

Note that it is not immediately evident that we can apply the transformation 
from Lemma[2]to {uk,Vk), since the transformation (u, v) = DE{r, b) is only defined 
for {r,b) G M and (r, 6) = DE~^{u,v) for {u,v) G respectively. Since we only 
know that (u,v) G we do not have uniform boundedness. However, u,v 

take values ±oo at most on a set of measure zero. Hence we know that (r, b) G A4 
a.e., which allows us to apply the variable transformation. 


We define E : M C ^MC (f,b) ^ (r,b) = DE-\u,v), 

where (u,v) is the unique solution in to the linear problem 

(24) 


a((u,'u), ($i,$ 2 )) = F($i,$ 2 ) for all ($i,$ 2 ) G 


with 


a((u,v), ($i,$ 2 )) = 


F($i,$2) = -- 


V$i 
o 

r - Tk-i 
b - bk-i 


G(r,b} 


Vu 

Vv 


$1 

4*2 


dx dy 


dx dy + tR 

- f H{P,b) 


$1 
$2 
V4>i 
V4>2 


dx dy. 


The bilinear form a : R^) x77^(H; R^) —>■ R and the functional E : id^(H, R^) —>■ 

R are bounded. Moreover, a is coercive since the positive semi-definiteness of G{r, b) 
implies that 


Vit 


a{{u, v), (m, v))= ( j G(f, &) ( j dx dy tR 


Vu 


> r M 


Iffi(a) 




Then the Lax-Milgram lemma guarantees the existence of a unique solution (u, v) G 
H^(H;R^) to ((Ml) . 

To apply Schauer’s fixed point theorem, we need to show that E is continuous. 
Therefore, let {rk,bk) be a sequence in M converging strongly to (f, b) in R^) 

and let {uk,Vk) be the corresponding unique solution to (IMll in 77^(H;R^). We 
have that G{fk,bk) —>■ G{f,b) and H{rk,bk) H{r,b) strongly in L^(fl,R^). As 
the entropy inequality yields a uniform bound for {uk,Vk) in i7^(H;R^), there ex¬ 
ists a subsequence with {uk,Vk) {u,v) weakly in H^(fl;R^). In order to identify 
(u, u) as the solution of (1241) with coefficients (f, 6), we first consider problem (IMll 
only for test functions in ($ 1 , 4 ) 2 ) G IT^’°°(H, R^). Here, the (weak) limit (u,u) is 
well defined. Then, the L°° bounds of G{fk,bk) allow us to consider the problem 
(l24l) for all ($i,$ 2 ) S i?^(H,R^) applying a density argument. So, the limit (m,u) 
as the solution of problem dH with coefficients (f, b) is well defined. 

In view of the compact embedding iJ^(H,R^) ^ L^(H,R^), we have a subsequence 
(not relabeled) with {uk,Vk) —>■ (u,u) strongly in L^(H,R^). Since the limit is 
unique, the whole sequence converges. Together with the property that the map 
from {u, v) to (r, b) is Lipschitz continuous (cf. Lemma[3]), we have continuity of E. 
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Furthermore, the compact embedding ^ L^(ri;]R^) gives the compact¬ 

ness of T. Combined with the property that T maps a convex, closed set onto 
itself, we can apply Schauder’s fixed point theorem, which assures the existence of 
a solution {r,b) G M to (l24)l with {r,b) replaced by (r, &). 

The convexity of E implies that E{(pi) — E{(p 2 ) < DE{(pi) ■ (tpi — Lp 2 ) for all 
e M. Choosing (pi = {rk,bk) and ip 2 = {rk-i,bk-i) and using DE{rk, bk) = 
{uk,Vk), we obtain 

■ (u!) dxdy > ^ JjEinM) - E{rk-i,bk-i)) dxdy. (25) 
Employing the test function ($i, < 1 ) 2 ) = {uk,Vk) in (1^ and applying (1^ . we obtain 

J^E{rMdxdy + rJ^ dr.^b,) dxdy 

+ r^i/(rfc,6fc) + < J^Ed-id-i) dx dy. 

(26) 


Using the entropy inequality resolving recursion (EHl) leads to 



(27) 


4.2. The limit r — ?> 0. Let {rk,bk) be a sequence of solutions to (1^51) . We define 
rr{x,y,t) = rkix,y) and br{x,y,t) = bk{x,y) for {x,y) G U and t G {{k - l)r, fcr]. 
Then {r^^br) solves the following problem, where (Tt- denotes a shift operator, i.e. 
{arrr){x, y, t) = rr{x, y.t-r) and [(Jrbr){x, y, t) = br{x, y, t - r) for r < t < T, 


1 

r 



frr- cr^rd 
\br — Ckdr) 



dx dy dt 



7o 


9a:?'r(l “ 
9x&r(l - 

dyXril 
dybr{l 


Pr') T XrdxPT\ 

Pt') T bxdxPr J 

Pr) T XxdyPx\ 
Pr) T bxdyPx J 


dxdydt 

■ i^dyi) 



(dy{rxbr){l - 

Ady{r.K)[l- 

(1 - Pr)rrVVr\ 
{l-pr)brVVb) 


Pr) H" '^T^T^yPr\ 
Pr) '^T^T^yPrJ 





dx dy + tR 




dt = 0, 
(28) 


for ($i(t), $ 2 (t)) G L^(0,T;i7^(U)). Inequality (l?7l) becomes 
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E{rr{T),br{T)) dxdy + Co [ [ {1 - Pt)\'^ dx dy dt 

Jo Jq 

Co f / (1 - Pr)|V\/67P + |Vv^l - Prp + \V Prl"^ dx dy dt 
Jo Jn 

dt< E{ro,bo) dx dy + TC. 

Jn 


(29) 


10 Jn 

.rj\ 

Jo 


The previous inequalities allow us to deduce the following Lemma. Note that 
from now on K denotes a generic constant independent of r. 

Lemma 4. (A priori estimates) There exists a eonstant K G R.’*', such that the 
following bounds hold: 


||-\/l — PrVi/7v||i2(0_T;L2(O)) + ||-\/l — Pt'^\/K\\L^( o,T-,L^(n)) < K, 

llVl ~ Pr\\L^{0,T;H^{n)) + W Pt\\{ 0,T-,H^ (n)) < K, (30) 
\/^(lbT||L2(o.T;ffi(n)) + ||'yT||L2(o_'r;Hi(n))) < K. 

The bounds in (l30l) together with the L°°-bounds for br and pr imply 


ll^^•^(l — Pt) + ||L2(o^T;L2(n)) 

<2||v^||L°°(0,T;L°°(n))||\/l - Pt\\ L“(0,T;L“(n))|| \/l — PrVv^||L2(o,T;L2(n)) 
+ lkT|U“(0,T;L“(n))||Vpr||L2(o,T;L2(n)) < K, 


(31) 


with an analogous inequality for br- Similarly, we get the estimate 

||V(r.r&r)(l ~ Pt) + l^rbr'^PT\\L^(0,T-,L^{n)) ^ (32) 

For applying Aubin’s lemma, we need one more property involving the time deriva¬ 
tives of Xr and br- 

Lemma 5. The discrete time derivatives of Xr and br are uniformly bounded, be. 

-\\rr — O'Trr||L 2 ( 0 ^T;ffi(n)') + “ll^-r — CTT&r||L2(o_T;i/i(n)') < K- (33) 
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1 

T 


Proof. Let $ G L^(0,T; Using the estimates in (I5(I1) . (1511) and (15^ yields 

■ [ {rr -crrrr,^)dt 
Jo 

■ / / [d^rril-P t)+ rrda:Pr)d^^dxdydt 

Jo Jn 

■T r 

(lodyTril - Pt) + rrdyPr)dy^ dx dydt 


h 

2 


— h 


Jo Jn 

/ {ldy{Trbr){^ — Pt) + T’TbTdyPr)dy^ dx dy dt 

Jo Jn 

i-T 


f ( (1 — Pt)?'tVU. • V<i> da; dy dt 

Jo Jn 

cT r 

Ut^ + Vut • dxdydt 


lo Jn 


— — Pt) + ^TdxPT\\L‘^{0,T-,L'^{n))\\dx^\\L^{0,T-,L'^{n)) 


+ ^||7o9yrT-(l — Pt) + ?"rl9yPT||L2(0_T;L2(n))||9y$||i2(0,T;L2(O)) 

+ b\\^dy{rTbT){l — Pt) + '<'TbTdyPT\\L^{o,T-L'^{n))\\dy^\\L'^{o,T-L'^(n)) 
+ 11(1 “ PT)?'rVK||L°°( 0 ,T;L“(n))||V<i)||il(o,T;Li(n)) 

+ ''■||wT||L2(0,T;_fai(n))||‘&||L2(o,T;ffi(n)) 

< l^ll‘J’llL2(0,T;Hi(n))- 


A similar estimate can be deduced for b which concludes the proof. 


□ 


From Lemma 0] and Lemma El we know that Pt G L'^{0,T\H^{Vt)) and ^{pt - 
(JtPt) G L^(0,T; respectively. This enables us to use Aubin’s lemma (cf. 

03], Theorem 1) to conclude the existence of a subsequence, also denoted by Pt, 
such that, as T —>■ 0: 

Pt ^ p strongly in L^(0,T;L^(r2)). 


This implies 

1 — Pr —t 1 — p strongly in L^(0,T;L^(r2)), (34) 

Pt ^ \/l- P strongly in L"^(0, T; L^(n)). 

Due to the continuous embedding of T; L^(n)) in L^(0, T ; L^(n)), it also holds 
that 

Vl - Pr -)• \/l - P strongly in L^(0,r;L^(r2)). (35) 

To pass to the limit t — 0 in we need to identify the weak L^- limiting func¬ 
tions of the following terms: 


(i) Vr7.(l — Pt) + Ci-Vpi-, V^t)! — Pr) + 6tVpt 


(ii) c^y(rTdT)(l Pt) F XTbTcJyPT 
(hi) (1 - pT)rTyVr, (1 - pT)bTyVb 
(iv) tUt, tVt, rVuT, tVut 






18 


Burger, Hittmeir, Ranetbauer, Wolfram 


The terms in dml) converge weakly in T^(0, T ; as l — pr converges strongly 

in L^(0, T; L^(n)) by ((Mil and because of the L°° bounds for br and r,-, up to a 
subsequence, 


Vr ^ r, br ^ b weakly* in L°°(0,T;L°°(f2)). (36) 

Because of O, we get that 

TUr , TVr 0 Strongly in (0, T; iJ^ (f2)), 
which identifies the limit in (liM . 

The weak convergence of the terms in ([il) and ([nl) can be shown with the help of a 
generalized Aubin-Lions Lemma (see Lemma 7 in [2D])- It states that if (IHHll . 

(l36l) and 

ll\/l - Pt 9\\L^{o,T-,H^n)) < K for g G {l,rr,br} (37) 

hold, we have strong convergence up to a subsequence for all f = f{rr,br) G of 

'/f-^f{rr,br) ^ \/T^f{r,h) strongly in L^(0,T; L^(n)), (38) 

as T —>■ 0. Note that (123 can be deduced from the previous a priori estimates. 
Writing (jl]) as 

Vr 7" ( 1 Pr ) T Tr ^Pt — Pt^ ( \/l Pt^t ) \/l Pt^ \/l Pr ^ ( 1 Pr 

— \/l Pt^ (\/l Pt^t ) \/l Pt^ \/l Pr 

- 2rr^/l - Pr^^/f - Pt 

= yjl- PrV(\/l - PrTr) - 3rr ^/l - PrVi/l - Pr 


and applying (|25)) with f(rr,br) = Tr, we get that 

\/l - PrTr ^ pr strongly in L^(0,r,L^(n)). 

Moreover, the L°° bounds together with (IHHI) give us bounds for V(i/1 — PtTt) = 
V-v/1 — — Pr'^y/fw and Vy^l — Pt- 

Together, we have 

VrT-(l — Pr) + rrV Pr Vr(l — p) + rVp weakly in L^(0,T; L^(r2)). 

Similarly, © can be written as 

^y(^r^r)(l Pr) T r" rbr^y Pr — Pr^yi^f Pr'^r^r) \/ f Pr^y \/ f Pr • 


Applying (I38p with f(rr,br) = Trbr and using analogous arguments as in ([I]), we 
get that 

dy{rrbr){f — Pr) + ^rbrSyPr dy{rb){f — p) + rbdyp weakly in L^(0,T; L^(ri)). 

(39) 

From Lemma [5] we derive that 

— {rr — (Jrrr)^dtr, —{br — <Jrbr)^dtb weakly in L^(0, T; 77^(0)'). 

T T 

This, together with the convergences in (I251) - (I2D1) allows to finally pass to the limit 
T —>■ 0 in (1^51) . which gives the weak formulation dlM) . 

The only thing which remains to verify is the entropy inequality (EDD- Since E is 
convex and continuous, it is weakly lower semi-continuous. Because of the weak 
convergence of {rr{t),br{t)), 


/ E{r{t),b{t))dxdy <\im inf / E{rr{t),br{t)) dx dy for a.e. < > 0. 

Jo. ^^0 Jo. 

We cannot expect the identification of the limit of \/l — pr'S/y/T\, but employing 
(1251) with f{r, b) = y/r, we get 


y/l- Pry/^ ^ 1 - pVr strongly in L^(0, T; L^(n)) 
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with analogous convergence results for r being replaced by b. Because of the L°°- 
bounds and the bounds in (23), we obtain V(-v/l — Pry/^) S L^(0, T; L^(n)), which 
implies 

\/l - ^ \/l- pVr weakly in 

— Pr\/K ^ v^l — weakly in L^(0,T;_ff^(ri)). 

The L°°-bounds, (HU and the fact that 

Va/1 — Pi- ^Va/1 — p weakly in L^(0,r;L^(n)), 


imply 

(1 Pi-) V-y/lv — \/l Pr^("x/1 Pr\/^T) \/l Pr'x/^r^^ \/l Pr 

(1 - P^)V\/^ = -v/l - PrV(v^l - Pr'/K) - y/l - PtVK'^V^ - Pr 

converge weakly in to the corresponding limits. The bounds imply also weak 
convergence in 

(1 — pi-)V-yrv ^ (1 — p)V-\/r weakly in L^(0, T; L^(n)), 

(1 — pi-)V'\/^ ^ (1 — p)vVb weakly in T^(0, T; L^(n)). 

As 1 — Pr > (1 “ Pt)^, we can pass to the limit inferior r —?> 0 in 





attaining the entropy inequality ([ 20 ]) • 

4.3. Existence for the general model. The previous analysis can easily be ex¬ 
tended to the case 0 < a < ^ and 71 — 72 = 0{h) in ([5]), i.e. leading to the system 
The particular choice of parameters allows to obtain the following result. 

Theorem 4. (Global existence) LetO < a < ^, 71—72 = 0{h), T > 0 and (ro,&o) : 
n ^ M, where A4 is defined by (US, be a measurable function such that E(rQ, bo) G 
L^(n). Then there exists a weak solution (r, 6) : fl x (0, T) —>■ At to system ® 
with periodic boundary conditions in x-direction and no-flux boundary conditions in 
y-direction satisfying the same regularity results and entropy dissipation inequality 
as stated in Theorem d 

The parameter regime 71 — 72 = 0{h) and 0 < a < ^ results in additional terms 
in the entropy dissipation which can be estimated using Young’s inequality in a 
way that inequality da holds for the same constant Cq. In the existence proof, 
the limit r —>■ 0 requires some additional compactness results for the new terms 
which we obtain (as in the proof of Theorem 13 by a generalized version of the 
Aubin-Lions lemma (cf. Lemma 7 in |30jl. 
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5. Numerical simulations 


In this last section we illustrate the behavior of the model with numerical sim¬ 
ulations in spatial dimension two. In particular we compare the solutions of the 
minimal model (HU, i.e. no cohesion and no preference for dodging to one side 
with those of the system ([S]). All simulations have been carried out using the 
COMSOL Multiphysics Package with quadratic finite elements. We consider the 
domain = [ 0 , 1 ] x [ 0 , 0 . 1 ] representing a corridor, where we use a mesh consisting 
of 608 triangular elements and a BDF method with maximum time step 0.1 to solve 
the corresponding system. We start with small perturbations of trivial stationary 
states to study if the system returns to this trivial solutions or results in a more 
complex one. 


5.1. Example I: Equilibration. We consider system m with no adhesion and 
no preference to step to the right or left, i.e the minimal model (HU). We choose 
the parameters 70 = 0.1, 7 = 0.2 and h = 0.3 and the initial values 


ro{x, y) = Cr + 0.02 sin 
bo{x,y) = Ch — 0.02 sin 


(Tra;) cos , 

(7ra:)cos , 


(41) 


with Cr = Cb = 0.4. Figure HI illustrates the initial value xq and the solution xt to 
system HU at time T = 5, where it can be seen that in this setting the solution 
returns back to the equilibrium state quickly. 


Tlme-0 s Surface; Dependent variable r (m) 



Tlme-5 s Surface: Dependent vanable r (m) 


0.1 0.2 0.3 0.4 


O.S 0.6 0.7 0.8 0.9 1 


(b) ry for T = 5 


A 0.4 


0.41 


0.4 


0.39 


▼ 0.4 


Figure 2. Example I: Red particle density returning to the con¬ 
stant stationary state after initial perturbation. 


5.2. Example II: Lane-formation. The behavior of solutions to system (lU cor¬ 
responding to the scaling 71 — 72 = 0{h) is different. Setting 70 = 0.001, 71 = 0.5, 
72 = 0.4, a = 0.2 and h = 0.1, such that 71 — 72 = 0{h), and choosing the same 
initial values (1411) as above, we obtain weak lane formation illustrated in Figure [3| 
As 7 i > 72 , the individuals have a tendency to step to the right. Therefore red 
individuals are highly concentrated on the bottom of the domain, whereas the blue 
individuals move to the top. 

Figure |4| shows the cross section of the two-dimensional solution xt at time 
T — 100 for different initial masses. The parameters are the same as above while 
the constants = C{, = c in the initial values m are being varied. We observe that 
weak lane formation is more pronounced for smaller values of 70 as well as higher 
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(a) vt for T = 5 


(b) bT for T = 5 


Figure 3. Example II: Red and blue particle distribution at time 
T = 5 forming weak lanes. 


densities. The transition region around y = 0.05 decreases for smaller diffusitivity 
and greater mass, while the behavior of r is the same in the low density region 
(0.05,0.1) for all parameter sets. 

In the case of different masses Mr and we observe asymmetric weak lane 



—c=0.25, Yo=0.0001| 


|_c=0.25, Yo=0.0001 



—c=0.25, Yo=0.001 1 


I— c=0.25, Yo=0.001 



—c=0.1, Y(,=0.0001 1 


|—c=0.1, Yq=0.0001 


X 

c=0.1, Y(,=0.001 


I c=0.1, Yq=0.001 






(a) r at time T = 100 (b) b at time T = 100 


Figure 4. Example II: Red and blue particle density in y-direction 
at time T = 100. 

formation. We choose initial values of the form gH), i.e. Cr = 0.4 and Cf, = 0.1. 
Figure[^shows the formation of such weak lanes due to the side-stepping mechanism 
even though the mass Mf, is smaller than Mr- 



Figure 5. Example II: Red and blue particle density as well as 
their sum p at time T = 100 in case of non-equal initial mass. 


5.3. Example III: Jam. We conclude with a numerical simulation showing an¬ 
other well known phenomena in crowd dynamics, namely traffic jams or so-called 
’freezing’. If the diffusion coefficients are small, i.e. /i = 0.1 and 70 = 0.0001, it 
may happen that the individuals cannot move in their walking direction any more 
as the initial masses are high compared to the diffusion coefficients. This ends in a 
jam or ’frozen’ configuration as Figure [5] illustrates. 
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(a) vt for T = 100 


(b) bT for T = 100 


Figure 6 . Example III: Congestion in the red and blue particle 
density resulting in a deadlock. 
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